arXiv: 1509.02318vl [cond-mat.stat-mech] 8 Sep 2015 


Stability of patterns on thick curved surfaces 

Sankaran N. * 

School of Physics, Indian Institute of Science Education and Research, Thiruvananthapuram-695016, India 

(Dated: September 9, 2015) 

We consider reaction-diffusion equations on a thick curved surface and obtain a set of effective 
R-D equation to O(e^), where e is the surface thickness. We observe that the R-D systems on these 
curved surfaces can have space- dependent reaction kinetics. Further, we use linear stability analysis 
to study the Schnakenberg model on spherical and cylindrical geometries. The dependence of steady 
state on the thickness is determined for both cases, and we find that a change in the thickness can 
stabilize the unstable patterns, and vice versa. The combined effect of thickness and curvature can 
play an important role in the rearrangement of spatial patterns on thick curved surfaces. 
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I. INTRODUCTION 

In 1952 Turing [T] proposed the reaction-diffusion 
(R-D) mechanism, where the chemicals can react and dif¬ 
fuse so as to produce spatially varying concentrations of 
chemicals in the steady state. Since then many mod¬ 
els HEl have been proposed to mimic the complex pat¬ 
tern formation in biological systems. Turing-like R-D 
equations are routinely used in trying to understand the 
skin patterns of animals [ 2 ] ; for example in fish | 3 ] , mam¬ 
mals [5], snakes [ 6 ] leopards [7] and many others. There 
have also been attempts to study changes in the pigmen¬ 
tation patterns on leopards and jaguars as they grow in 
size [7]. 

In most studies, R-D equations are analyzed on flat 
geometries which are not always suitable for the study of 
patterns on animal skin surfaces. It is reasonable to as¬ 
sume that geometry of the surface can play an important 
role in determining the pattern formation. For instance, 
Turing considered the surface of a sphere in the context 
of gastrulation of a blastula [T]. Geometry is probably 
responsible for stripes at the end of the tail while there 
are spots elsewhere in some animals [2]. Understand¬ 
ing the pattern formation on curved surfaces can be im¬ 
portant in some chemical, biochemical and embryological 
process [5]. It is also known that organ morphogenesis 
can be controlled by tissue geometry [3]. Recently some 
studies have been initiated to understand the role of cur¬ 
vature in biological systems nnHisi- 

The thickness of the surface is being neglected in 
most earlier studies. For instance, blastula is considered 
as a hollow sphere assuming no thickness [T] . But protein 
diffusion in lipid bilayers can be viewed as a diffusion on a 
two-dimensional curved surface with thickness. In some 
of the recent studies, there has been attempts to incor¬ 
porate the small thickness [Tiins] . The combined role of 
geometry and thickness can lead to curvature-dependent 
diffusion and may result in complex pattern formation 
in animals M- Recent study on the pattern formation 
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in melanocytic tumours again suggests the importance of 
geometry and of the thickness [15]. 

Some of these studies suggest that it is relevant to 
ask about the effect of the thickness and curvature in R-D 
systems. In particular, how does the curvature and thick¬ 
ness affect the formation of steady state patterns?. We 
answer this question by considering the reaction-diffusion 
of two chemicals on a curved surface with small thickness. 
We first obtain an effective description of R-D equation 
and then deduce the dependence of steady state and its 
stability on the thickness. 

This paper is organized as follows. In sec. II, we 
describe a general model of a R-D system and then ex¬ 
plicitly obtained its effective description. In sec. Ill, 
we analyze the effect of the thickness and curvature in 
the Schnakenberg model, specifically on a spherical and 
cylindrical geometry. Finally, we summarize our results 
in sec. IV. 


II. MODEL 

We consider reaction-diffusion of two chemicals be¬ 
tween two curved surfaces, a and cr , which are parallel 
to each other and separated by a distance e. The concen¬ 
trations of chemicals are denoted by A{q°, , t) and 

B{q^,q^,q^,t) and the dynamics is governed by the R-D 
equations [T] 


Fi{A,B) + DaV^A 

( 1 ) 

F 2 (A,R) 

( 2 ) 


where ^ 1 (^ 1 , B) and ^ 2 ( 2 !, B) are the the reaction kinet¬ 
ics, which in general are nonlinear functions; Da and Db 
are the diffusion constants of the chemicals. 

The co-ordinate system we use here is similar to the 
one considered by Ogawa M- We place a curved surface 
U between the surfaces a and tr at distances e /2 and 
—e/2 respectively. Any point between the surfaces a and 
a can be represented by q^ = {q^,q^,q^) where {q^,q'^) 
are the curved co-ordinates on the surface D and < 7 ° is the 
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FIG. 1: Coordinate System 


normal co-ordinate. The components of the metric 
of the curvilinear co-ordinate system can be given as 


Goi — Gio — 0, Goo — Ij 

Gij = dij + Kij + KimKj^ , (3) 

where gij is the metric on the curved surface 17, Hij is 
the second fundamental tensor. The determinant of the 
metric, G = det(G^y) can be written as, 

G =g {l + 2q°K+iq^f{^ + R)+0{e^)}, (4) 


where 


and 


FiiAB) = 
F2iAB) = 


A = 


eff 




.e/2 

-e/2 

j.e/2 

/-e/2 

ee/2 

J-e/2 

.e/2 

'-e/2 


dq\ - FiA,B), 


dq\ - G{A,B), 


dq° 


dq° 




g 

[g 




( 11 ) 

( 12 ) 

(13) 

(14) 


If we assume that concentrations of chemicals are inde¬ 
pendent of co-ordinate, then we obtain to 0{e^) 


A = eil+—R) A{q\q‘^), 
B = e{l + ^R) B{q\q^), 


(15) 

(16) 


and hence equations (11) and (12) can be rewritten as 


where, g = det(gij), mean curvature k = g^^Kij and Ricci 
scalar R = {k^ — 


A. Effective Theory 


In this subsection, we obtain the effective two di¬ 
mensional description of R-D equations. 

The total amount of the chemicals present in the 
system can be decomposed as follows 


J AVg d^q = 

[ bVG d^q = 



Vg d?q, 

Vg 


(5) 

( 6 ) 


thus leading to a definition of concentrations of chemicals 


in the effective description 



.e/2 

to 

II 

/ dq^ 


l-e/2 


.el2 

B{q^,q^,t) = 

/ dq^ 


l-e/2 


^ A{q°,q\q'^,t), (7) 

B{q°,q^,q^,t), ( 8 ) 


Fi(i,R) = e(l+^i?)Fi(J(l-|^R)i,i(l-|^i?)R), 

(17) 

F2(i,B)=e(l + ^i?)T®2(J(l-|^i?)i,^(l-^i?)S). 

(18) 

Assuming the fluxes VA and VR vanish at boundaries 
a and cr , and following similar steps as in [14j for the 
equations (13) and (14) will lead to 


Vlff = a ( 2 )+v 


(19) 


V = g 

12 


1/2 


dq^ 


X {(3/1 


d 


1 


- 2 /^/^®^)—}. 


dq^ 


dqi 


where is the Laplace-Beltrami operator on the 

curved surface 17. To summarize, the effective descrip¬ 
tion is captured in (9),(10),(17),(18) and (19). 


III. EFFECT OF THICKNESS AND 
CURVATURE 


Multiplying equations (1) and (2) with \fGjg and inte¬ 
grating over result in the equations 


h{A,B) + DA^liiA, 

(9) 

i?2(i,R) + DBV2^^R, 

(10) 


In this section, we illustrate the effect of the thick¬ 
ness and curvature in the effective description of a R-D 
equation. In particular, how does the nature of steady 
state vary by changing the thickness?. This question can 
be addressed within the framework of above discussed ef¬ 
fective theory. Note that, the thickness is kept constant 
during the dynamics of the system. 
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In flat geometries, it can be seen from equations (15) 
and (16), both the concentration are scaled, A eA and 
B —>■ eB, and dynamics remains the same. While in 
curved geometries, combined role of thickness and curva¬ 
ture can lead to nontrivial effects. 

In the following, we consider a R-D system between 
two curved surfaces, and point out some salient features 
of the effective description. We then proceed to analyze 
the effective description on spherical, and cylindrical ge¬ 
ometry. 


A. General surface 

Let us restrict to the well-studied Schnakenberg model, 
where the reaction kinetics is taken as [16] 

Fi{A,B) = ki-k2A + k3A^B, (20) 

F2{A,B) = ki-k3A^B, (21) 

Here, we confine the chemicals between two curved 
surfaces a and a placed at a distance | from a general 
curved surface H. We can read the effective description of 
Schnakenberg model on a general surface using equations 
(17) and (18), and is given by 

r)A 

= ki — k2A k^A^B + DA^'lffA, (22) 
dB 

— = h-hA^B + DBVljfB, (23) 

where given by (19), and the reaction constants 

in effective theory are related to the original model as 

ki = €ki{l + ^R), k3 = ^ksil- ^R), 

£2 

ki = eki{l + k2 = k2 ■ 

In general, the reaction constants in an effective the¬ 
ory of the Schnakenberg model are space-dependent as 
the Ricci scalar(i?) is not necessarily a constant, and may 
lead to interesting consequences. The space-dependent 
reaction kinetics can result in an absence of a homoge¬ 
neous steady state. 

A few comments about the dependence of reaction 
rates on Ricci scalar follows. The term k^A^B in equa¬ 
tion (22) represents the production of the chemical A. 
Note that the reaction constant k^ is lower in regions 
with higher positive curvature. Hence, the production of 
A is more in regions with lower positive curvature. Sim¬ 
ilarly the term -k^A^B in equation (23) represents the 
depletion of B and hence can result in more depletion 
in regions of lower positive curvature. Both the reaction 
constants ki and ki are higher in regions with higher posi¬ 
tive curvature and result in more production of chemicals 
in these regions. 


The term DAS/^ff can be rewritten as 

-^d,y^{Dp,A - ^^g^^d.R.A), (24) 


where 


DJ = DA{g^^ + - 2kk^^}). (25) 

The first term in equation (24) is the diffusion term, 
which is not necessarily isotropic, and is characterized 
by the diffusion matrix which depends on the extrin¬ 
sic curvature. On a sphere the O(e^) term of both 
and is negative, where (0, (p) are the co-ordinates on 
the surface of a sphere. But on a cylindrical surface the 
0{€^) term of is positive and = 0, where {9,z) 
are the co-ordinates on the surface of a cylinder. Hence 
there is an enhanced diffusion of chemicals along the the 
9 direction on a locally cylindrical region. The diffusion 
of chemical can be slow down along 9 and (j) directions on 
locally spherical regions. The second term in the equa¬ 
tion (24) is the current due to the gradient of Ricci scalar 
between two points on a surface [T3]. 


B. Spherical Geometry 

We now consider the Schnakenberg model, where the 
chemicals are conhned to the region between two spheres 
of radii Oq -f | and Oq — |. Choose the sphere with radius 
oo as the surface and spheres with radii oq -b f and 
oo — § are the surfaces a and cr, respectively. Our 
analysis suggests the following effective description for 
this model 


dA 

~di 


ki — feA + it^A?B -b VgjjA, 
f) n 

— = fc4 - hA^B + ^IffB, 


(26) 

(27) 


where the reaction constants in the effective theory are 
related to original model as follows 


ki = kie{\ -b 


24al^’ ^ 2 -^ 2 , 


since R — I /og for the sphere of radius og. 

Choosing the {9, (j)) coordinates on the surface of a sphere 
with radius og, it is straightforward to obtain 


gee = Og, 300 = Og sin^ 9, gg^ = g^g = 0, 

f^QQ — ^0; — ^0 sill ^9(f> — ^cj)9 — 0; 

4> ^ 9 ^ <p 9 r\ 

Kl = —= = 0 . 

ao ao 
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Hence read as 


The solution to the equation (32) can be written as 


2 _ 1 92 

~ 62 ^ 902 




where bo = ao(l + 23 ^)- essence, the effective two 
dimensional R-D equations (26) and (27) can be inter¬ 
preted as Schnakenberg model on a sphere of radius bo 
with redefined parameters. 

Equations (26) and (27) can be rewritten in terms of 
rescaled variables as 


r)TI 

_ = ^(«_C/ + y2 V) + ^lffU, 
8V 

V) + dV%fV, 


(29) 

(30) 


where. 


r = DAt/bl U = i(^)l/^ E = H(^)i/2, d=^, 
k2 K 2 


(31) 


/C2 ^2 k,2 f^2 


blK2 


and is the Laplace operator in the scaled variables. 

The homogeneous steady state solution can be ob¬ 
tained from (29) and (30) as {Uo,Vo) = {a + b, ). 
We consider only positive solution of homogeneous steady 
state and discard the solution with negative concentra¬ 
tion. Note that homogeneous steady state solution is 
independent of the thickness to C>(e2). 

The linear stability analysis about the homogeneous 
steady state follows. A small variation in the homoge¬ 
neous steady state is denoted as 


5W = 


dU-Uo\ 
6V-Vo )' 


which satisfies the linearized equation 


where 


d {5W) 
dt 


= LSW, 


L = jC + DVlff, 


(32) 


(33) 


D = 




dU dV 
dg dg 

au dV 


Uo,Vo 


and 


/(U, E) = 7(a - 17 + C/2 E), 
5(C/,E)=7(C-C/2 E) . 


l—oo I 

SW{e,(j),t) = Y^ ^ C're^(')*P/"(cos6>)e*’"‘^, (34) 

l—O rn— — l 


where the constants (7[" can be determined from initial 
conditions. The eigenvalues X{1) satisfy 

+ X[{l{l + l)){l + d) — 'y(fu+gv)] + h{l{l + l)) = 0, (35) 
where fu = §ij,gv = -§y and h{l{l + 1)) can be given as 
h{l{l+l)) = d(l{l+l))^—'j{d fjj_+gy) l{l+l)+j'^{fugv—fvgu) 

(36) 

The necessary condition for the instability to kick in is 
h{l{l + l)) < 0, 


The modes I which satisfy h{l{l -|- 1)) < 0 gives 
positive eigenvalues in equation (35). These are the 
modes which give rise to the instability (an inhomoge¬ 
neous steady state) to the system. These modes (unstable 
modes) lie between L_ < 1(1 -|- 1) < L+, where L_ and 
L+ are the roots of li(l(l -I- 1)) = 0, and given as 



FIG. 2: h(x) Vs x: Figure illustrates the role of thickness- 
dependent 7 in the selection of unstable modes. 


Note that the unstable modes lie between the 
zeros(L+, L_) of the function h{x), where L+ and L_ de¬ 
pend on the thickness. It can be seen from the equation 
(37) that an increase (decrease) in thickness can shift the 
zeros of h{x) towards right (left) on the x axis as shown 
in the fig(2), and this can result in a scenario where the 
stable modes can become unstable and vice versa. Hence 
the nature of steady state can be changed by tuning the 
thickness. In other words, there can be a transition from 










5 


a homogeneous steady state to an inhomogeneous steady 
state(a pattern) by changing the thickness. It is also pos¬ 
sible to obtain different inhomogeneous steady states as 
thickness changes. The e dependence in the equation (37) 
thus reveals the effect of thickness on the stability. 

Let us distinguish two cases. In the case I, for certain 
ranges of parameters(a,b,d, 7 ) the zeros of function h{x) 
lie below the point x = 2 , namely, 0 < L_ and L+ < 2 . 
Here, the steady state configuration is homogeneous. In 
this case, since there is no unstable mode lies below the 
point X = 2 , the instability can set in only by increasing 
the thickness. 

In the case II, the zeros of function h{x) lie between 
li{li + 1) and (^1 -I- l)(Zi -I- 2), namely, li{li -I- 1) < L- 
and L+ < {li + l)(Zi -I- 2). In this case the system can 
be driven to an inhomogeneous steady state either by 
decreasing or increasing the thickness. 


C. Cylindrical Geometry 


where the scaled variables are defined in equation (31) 
and z = and is Laplace operator in the scaled 

variables. 

If we now proceed similar to the case of a sphere, 
then the deviation from homogeneous solution 

SW{e,z,t) = ^C'„me^*e™®cos("^*^z), (43) 

n,m 


where Cnm depends on the initial conditions. The eigen¬ 
values can be obtained from 


X^ + X[{ri^ + - 


.2^2.2 


P 


){l + d)-'y{fu+gv)]+h{n,m) =0, 


where 


h{n, m) = d{n'^ -\- 


2 


P 


} - l{d fu+gv){n + 


m TT 


P 




T7 (/itffr fvgu)- 


Here we consider Schnakenberg model, where the 
chemicals are confined between two cylinders of radii 
oo -f I and ao — |, and further assume the flux vanishes 
at z = 0 and z = 1. 

In this case the effective description is governed by 

dA _ . _ 

= Ki — K2A-\-fA, (38) 

f) Pi 

— = ki-ksA^B + VlffB, (39) 

where the reaction constants in the effective theory are 
related to the original model as follows 

ki = ki€, k^ = ^k^, k4 = kit, ^2 = ^2 ■ 

on the surface of a cylinder with radius oq the quantities 
related to intrinsic and extrinsic curvatures are 


gee = ai,gzz = = gze = 0, 

i^ee = 0 , 0 , = 0 , Kgz = Kzb = 0 , 

= 0 . 

ag 

Hence, 

2 


(40) 


where 60 = ao{l ~ Thus the effective equations 

(38) and (39) can be interpreted as Schnakenberg model 
on a cylinder with rescaled radius &o- Equations (38) 
and (39) can be rewritten in terms of rescaled variables, 
U = U(0,z), and V = V{9,z), obeying the equations 


_ = ^(«_C/ + y2 V) + VlffU, 
V) + dVlffV , 


(41) 

(42) 


Following the same analysis as done in the case of a 
sphere, the modes (n,m) which satisfy the following con¬ 
dition 


(44) 

will destabilize the homogeneous solution, where BA and 
are the zeros of h(n,m) and given by 


Rlie) 





)[{dfuFgv) 


±{{dfu + gvf - 4:d{fugv - fvgAV^'^] ■ 


(45) 


Hence the unstable modes (n, m) lie between the two 
semicircles in the (n, —^—^) plane with radius R+ and 


R- as shown in the figure 3. It can be seen from equation 
(45) that the values of both R_ and R^. decrease by 
increasing the thickness(e). 



FIG. 3: Unstable modes (n, m) lie between semicircles with 
radius R- and R+. 
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Let us distinguish two cases. In the case I, the space 
between two semicircles in the (n, ) plane encloses 

no mode (n, m) satisfying the condition (44). The steady 
state configuration of R-D system is homogeneous in this 
case. By increasing the thickness the semicircles can 
shrink, and can result in the inclusion of some modes 
(n, m). Hence, there is a possibility of obtaining an inho¬ 
mogeneous steady state from the homogeneous one by in¬ 
creasing the thickness. It is also possible in this case that 
the region between two circles can encompass no modes 
(n, m) even after increasing the thickness. In such a situ¬ 
ation, the system can continue to be in the homogeneous 
steady state. 

In the case II, some modes (n, m) are enclosed within 
the semicircles of radius R+ and i?_. The steady state 
configuration is inhomogeneous in this case. Now an in¬ 
crease in thickness can result in the shifting of curves such 
that modes (n, m) are no longer present between i?_ and 
i?+ semicircles. Hence, there is a possibility of transition 
to a homogeneous steady state from an inhomogeneous 
one by increasing the thickness. Another possibility is 
that the semicircles can shrink. This can lead to an in¬ 
clusion of new modes (n, m) and result in transition to 
a new inhomogeneous steady state from the initial inho¬ 
mogeneous state. 

IV. CONCLUSION 

To conclude, we have studied the effect of curva¬ 
ture and thickness in R-D systems on quasi-two dimen- 
sional(thick) curved surfaces. We explicitly analyzed an 
effective description of the Schnakenberg model, in par¬ 
ticular, on spherical and cylindrical geometry. In both 
spherical and cylindrical case, the effective theory is same 
as the original model on a sphere and a cylinder, respec¬ 
tively, with rescaled parameters. On the spherical geom¬ 
etry an increase in the thickness can lead to an increase in 
the parameter 7 . In cylindrical geometry the parameter 
7 can decrease by increasing the thickness. In the ab¬ 
sence of curvature, the thickness play no significant role 
in the effective description. 

In general, R-D systems on quasi-two dimensional 
curved surfaces can have space-dependent parameters. 
There are a few instances where spatially varying pa¬ 
rameters are considered The absence of homo¬ 

geneous steady state is also a characteristic of the ef¬ 


fective R-D equation. The effective R-D description is 
not necessarily similar to the original description on a 
two-dimensional surface with rescaled parameters. In 
R-D systems on quasi-two dimensional curved surfaces, 
a change in thickness can stabilize the unstable pat¬ 
terns and vice versa. Hence the patterns (inhomoge¬ 
neous steady state) can appear or disappear by tuning 
the thickness. This might be a plausible reason for dif¬ 
ferent patterns on leopards and jaguars as they grow 
in size or rather as the skin thickness increases. 

There is a related model studied in the context of 
Belousov-Zhabotinsky reaction [52] where there is no dif¬ 
fusion. Instead, equation to the chemical A contain v.VA 
term, where v is the velocity of the chemical A, while the 
other chemical B is immobilized. In this model the chem¬ 
ical instability is of traveling-wave type, and the concen¬ 
trations can vary both in space and time. In this case 
there is an isotropy in the wave speed when the speed 
of the chemical A is same in every directions. Assume 
that the velocity(z;) of the chemical A is independent of 
the < 7 ° direction. Then following the methods described 
in sec.H, one can straightforwardly obtain the effective 
description of the above model on an infinite cylinder. 
The stability analysis follows provided the linear terms 
in the reaction kinetics of the effective theory meet the 
stability condition. The above analysis shows that the 
thickness can induce an anisotropy in the speed of the 
traveling-wave. 

The effective description outline in the paper can be 
easily extended to any R-D models like Gierer-Meinhardt 
model [23] , and other R-D models [21 El El] ■ The analy¬ 
sis may prove useful in the study of the rearrangement of 
spatial patterns during various stages of growth in ani¬ 
mals. Turing-like models also find applications in diverse 
areas like material sciences [23], hydrodynamics [23], as¬ 
trophysics [22], etc. In these systems, under certain con¬ 
ditions it is conceivable that the thickness and curvature 
can play a significant role, and hence similar effective 
descriptions may be suitable. 
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